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Abstract 

A Jfe — € model is proposed for wall bounded turbu- 
lent flows. In this model, the eddy viscosity is charac- 
terized by a turbulent velocity scale and a turbulent 
time scale. The time scale is bounded from below 
by the Kolmogorov time scale. The dissipation rate 
equation is reformulated using this time scale and 
no singularity exists at the wall. A new parameter 
R = ^ is introduced to characterize the damping 
function in the eddy viscosity. This parameter is de- 
termined by local properties of both the mean and 
the turbulent flow fields and is free from any geome- 
try parameter. The proposed model is then Galilean 
and tensorial invariant. The model constants used are 
the same as in the high Reynolds number Standard 
Jfe-e Model. Thus, the proposed model will be also 
suitable for flows far from the wall. Turbulent chan- 
nel flows and turbulent boundary layer flows with and 
without pressure gradients are calculated. Compar- 
isons with the data from direct numerical simulations 
and experiments show that the model predictions are 
excellent for turbulent channel flows and turbulent 
boundary layers with favorable pressure gradients, 
good for turbulent boundary layers with zero pres- 
sure gradients, and fair for turbulent boundary layer 
with adverse pressure gradients. 

1. Introduction 

In turbulence modeling, the fc-e model is the most 
widely used model in engineering calculations. The 
Standard k — e Model 1 * 2 was devised for high Reynolds 
number turbulent flows and is traditionally used in 
conjunction with wall functions when it is applied to 
wall bounded turbulent flows. However universal wall 
functions do not exist in complex flows and it is thus 
necessary to develop a form of k — e model equations 
which can be integrated down to the wall. 

Jones and Launder 3 were the first to propose a low 
Reynolds number k — e model for near wall turbu- 
lence, which was then followed by a number of simi- 
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lar k — e models. A critical evaluation of the pre-1985 
models was made by Patel et al. 4 . More recently 
proposed models axe found in Shih 5 and Lang and 
Shih 6 . Three major deficiencies can be pointed out 
about existing k — e models. (Some of the models 
may have only one or two of the three deficiencies.) 
First, a near wall pseudo-dissipation rate was intro- 
duced to remove the singularity in the dissipation rate 
equation at the wall. The definition of the near wall 
pseudo-dissipation rate was quite arbitrary. Second, 
the model constants were different from those of the 
Standard k — e Model, making the near wall models 
less capable of handling flows containing both high 
Reynolds number turbulence and near wall turbu- 
lence, which is often the case for a real flow situation. 
Patel et al. 4 put as the first criterion the ability of 
the near wall models to be able to predict turbulent 
free shear flows. Third, the variable y + is used in the 
damping function of the eddy viscosity formulae. 
Since the definition of y + involves u T , the friction ve- 
locity, any model containing y + can not be used in 
flows with separation. In addition, y + may not be 
well defined for flows with complex geometry. 

In an earlier paper by the authors 7 (referred as YS 
here after), a time scale based k — e model for near 
wall turbulence was proposed. In this model, k 1 / 2 
was chosen as the turbulent velocity scale. The time 
scale was bounded from below by the Kolmogorov 
time scale. When this time scale is used to reformu- 
late the dissipation rate equation, there is no singu- 
larity at the wall and the introduction of a pseudo- 
dissipation rate is avoided. The model constants were 
exactly the same as those in the Standard k— e Model, 
which ensures the correct performance of the model 
far from the wall. The damping function in YS was 
proposed as a function of Ry = - ! y instead of t/ + 
Thus, the model can be used for flows with separation 
and reattachment. 

However, the Ry dependence in the damping func- 
tion makes the model coordinate dependent. It also 
creates some ambiguity when the model is used for 
complex geometries, for example, a corner flow. Sim- 
ilar problem also exists in most existing k — e models. 



The aim of the present paper is to remove this de- expression for turbulent time scale, k/e, approaches 
ficiency while keeping the good performance of YS. to zero due to the boundary condition for k. How- 
This is achieved by introducing a new parameter R, eve i, as it was shown in YS based on a local Taylor 
to be defined and discussed in section 2, as the in- series expansion of the instantaneous velocity field, 
dependent variable of the damping function in the the turbulent time scale should approach to the Kol- 
eddy viscosity. The introduction of R also makes mogorov time scale near the wall rather than zero, 
the proposed model Galilean and tensorial invariant. The turbulent time scale given by equation (5) has 
The proposed model is then calculated for turbulent the properties that the turbulent time scale is given 
channel flows at different Reynolds numbers and for by the Kolmogorov time scale near the wall and the 
turbulent boundary layer flows at zero pressure gra- turbulent time scale is given by the conventionally 
dient, favorable pressure gradient, and adverse pres- used expression {k/e) far away from the wall, where 
sure gradient. The numerical aspects of the compu- Kolmogorov time scale is much smaller than k/e. 
tation are briefly discussed in section 3. The results k and e are found from solutions of the following 
of the model calculations and the comparisons with modeled transport equations 


the data from experiments and direct numerical sim- 
ulations are shown in section 4. Section 5 concludes 
the paper. 

2. Near wall k - e model 

In turbulence modeling, the instantaneous quanti- 
ties of an incompressible flow are decomposed into 
the mean and the fluctuating parts, i.e., Hi = Vi +Uj, 
p — p +p. The mean field Vi satisfies the following 
continuity equation and Reynolds averaged Navier- 
Stokes equation: 

Ui,i = 0 ( 1 ) 

Vi + Vj Vij = - 1 -. Pi + vVi'jj - u Wij ( 2 ) 

where the Reynolds stress term, u^uj) must he mod- 
eled. 

In an eddy viscosity model, one assumes that the 
Reynolds stress is related to the mean field hy 

— UiUj = vr{Ui,j + Uj t i) — 00 

where x/j * is the eddy viscosity and k is the turbulent 
kinetic energy* 

In the near wall k - e model proposed by YS, the 
eddy viscosity is written as 

v T = CpfpkTt ( 4 ) 

where 

T t = * + (^) 1/2 (5) 

is the time scale for turbulent flows. The first term 
on the right hand side of the above equation is the 
conventionally used expression for the turbulent time 
scale. The second term is the Kolmogorov time scale. 
As the wall is approached, the conventionally used 


El = [(»>+ — )M,j - wi u i,j - e ( fi ) 

Dt <Tk 

De rr . \ i ■ [-CuuiujVi'j - C 2f e) 

- = [(■'+-)..,■]„ + jr 

+ C t vi/'rVi i jhUi l jk‘ ( 7 ) 

The turbulent time scale given by equation (5) is 
used to reformulate the dissipation rate equation (7). 
Since the turbulent time scale is always positive, the 
dissipation rate equation does not have any singu- 
larity as the wall is approached* The argument on 
the time scale for near wall turbulence and the use of 
the Kolmogorov time scale for near wall turbulence 
modeling are also discussed by other researchers, see 
Durbin 8 , for example. 

f in equation (4) is the damping function which 
is introduced to account for the effect of the wall on 
the eddy viscosity. In YS, / M was given as a function 
of Ry - k^ 2 y/v. The parameter Ry has the advan- 
tage over the commonly used = u T ylv in that it 
could be used for flows with separation and reattach- 
ment* However, the appearance of y in the damping 
function makes the model non-tensorial invariant. In 
addition, in complicated geometry situations such as 
corner flows, the meaning of y is ambiguous. 

To overcome this deficiency, we introduce a new 
parameter R, which is defined as 

(8 > 

where S is the modulus of the strain rate tensor Sij 
of the mean velocity field and is given by 

S = (2 SijSij) 112 . ( 9 ) 

The parameter R defined above is expressed in 
terms of the local variables of the field and is thus 
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Figure 1 1 Variation of R for turbulent channel flow 
at Re T - 180. 


coordinate independent. The physical meaning of R 
can be explained as follows, 


k k 2 .Sk 

r =Tv = tJt' 


( 10 ) 


R is then the ratio of two important physical param- 
eters: the turbulent Reynolds number and the time 
scale ratio of the turbulence to the mean flow. The 
variation of R with the wall distance y+ for the case 
of turbulent channel flow at Re r = 180 is shown in 
figure 1. The gradual and monotonic increase of R 
with y+ in the near wall region makes jR an ideal 
candidate for constructing the damping function. 

In the present study, the form of the damping func- 
tion is chosen as 


which gives /„ -♦ 0 as O(y). Thus, the near wall 
asymptotic behavior for the shear stress is satisfied. 
Far from the wall, R is large and / M — ► 1. 

The other effect in the modeling of the near wall 
turbulence is the effect of the inhomogeneity of the 
mean field which introduces a secondary source term 
in the dissipation rate equation. This is represented 
by the last term on the right hand side of equa- 
tion (7). A term of this form in the dissipation rate 
equation was suggested by Jones and Launder 3 , and 
Shih 5 . The effect of this term is confined to the near 
wall region and is most prominent in the buffer layer. 
Away from the wall, this term becomes much smaller 
than the other terms in the dissipation rate equation. 

Equations (4), (6), and (7) along with the eddy 
viscosity given by equation (11) are the k — e equa- 
tions proposed in the present paper. The model could 
be used for complex flows, flows with separation and 
reattachment, for example. Since all the quantities 
are expressed in the local variables, the model can be 
incorporated in a general purpose CFD code, partic- 
ularly a code with unstructured grid. 

The model constants, c M , Ci e , C 2 « } <7*, <x 6 , are cho- 
sen to be the same as those in the Stanford k — e 
Model, i.e. — 0.09, Cu = 1*44, C 2 e = 1.92, 
a h = 1.0, <r € = 1.3. The value of C € is chosen to 
be equal to 1.0 based on the model performance for 
turbulent channel flows. Away from the wall, the 
present model reduces to the Standard k — € Model. 
Thus, it is only necessary to assess the performance 
of the model for wall bounded flows. 

The boundary condition for € on the wall is deter- 
mined by applying equation (6) at the wall, which 
gives 


— l/k^yy 

In this study, the following boundary condition for e, 
which is mathematically equivalent to the above but 
computationally much more robust, is used. 


fp — [l — exp(-aiR - a 2 R 2 — asR 3 )] 1 ^ 2 (11) 

where a x = 3xlCT 4 , a 3 = 6xl0“ 5 , a 5 = 2xl0“ 6 . 
These constants are obtained by calibrating the 
model predictions for the turbulent channel flow at 
Re r = 180 against the direct numerical simulation 
date of Kim et al. e 

Near the wall, the shear stress —uv should behave 
as 0(t/ 3 ), according to the local Taylor series analy- 
sis. Since k is 0(y 2 ) and the time scale T t is finite, 
we would require the damping function to have a near 
wall behavior of O(y). From equation (8) and equa- 
tion (11), it is seen that as y -► 0, R — ► 0 as 0(y 2 ), 


< 12 > 

3. Numerical aspects 

Boundary layer approximation is used in the cal- 
culations shown below. An implicit finite difference 
scheme is used to solve the momentum equation and 
the transport equations for k and e. The coefficients 
for the convective terms are lagged one step in the 
marching direction and the source terms in the k and 
e equations are linearized in such a way that numer- 
ical stability is ensured. 
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Figure 2: Computation of turbulent channel flow at 
Re r = 395 with different number of grid points. 

A variable grid spacing is used to resolve the sharp 
gradient near the wall. The grid distribution is con- 
trolled by 8yi/&yi-\ = a. Both a and the total num- 
ber of the grid, N } are varied to ensure the grid in- 
dependence of the numerical results. The marching 
step size, <5x, is also varied to ensure accuracy. In the 
calculations shown in the next section, these parame- 
ters axe changed such that the solution has a less than 
1.0% error. Typically, the grids used are specified by 
N = 150 and a = 1.05. It is also found that the solu- 
tion is not sensitive to the number of the grid points 
as long as there are two points in y+ < 1. Figure 
2 shows calculations for the turbulent channel flow 
at Re r = 395 with N varying from 50 to 150 and a 
changes accordingly such that a N remains basically 
the same. It is seen that the results for different N 
are almost identical. 

4. Results and discussions 

Turbulent channel flows at different Reynolds num- 
bers and turbulent boundary layers with zero pres- 
sure gradient, favorable pressure gradient, and ad- 
verse pressure gradient axe calculated using the 
present model. The following shows the computa- 
tional results along with the available experimen- 
tal data and data from direct numerical simulations. 
For some cases, the predictions of the Jones-Launder 
model and the Chien’s model^ are also shown. These 


two models are chosen because the Jones-Launder 
model is the first k - € model for near wall turbu- 
lence while the Chien’s model is known to perform 
quite well for turbulent boundary layer flows. 

Two dimensional fully developed channel flows 
were calculated first. These flows are attractive for 
model testing because they have self-similar solutions 
so that the initial conditions do not have to be accu- 
rately specified. These flows are very simple and solu- 
tions can be found very efficiently; yet, the effects of 
the wall on turbulent shear flow axe still present. In 
addition, DNS data providing detailed flow informa- 
tion are available for comparisons. Computations are 
carried out for 2D fully developed turbulent channel 
flows at Re r = 180 and Re T = 395, respectively. Fig- 
ures 3-6 show the profiles of the mean velocity, shear 
stress, turbulent kinetic energy , and the dissipation 
rate, respectively, for the case of Re r = 180. Both the 
dependent variable and the independent variables are 
represented in wall units by normalization through u T 
and i/. The predictions of the Jones-Launder model 
and the Chien’s model are also shown. These predic- 
tions are compared with the DNS data. It is found 
that the present model gives the right location of the 
maximum value of the dissipation rate. The center 
line velocity, which corresponds to the skin friction 
coefficient, is also well predicted. The corresponding 
results for the case of Re r = 395 are shown in figures 
7-10. 

Like the 2D fully developed channel flows, zero 
pressure gradient (ZPG) turbulent boundary layer 
flows over a flat plate also give a self-similar solu- 
tion. Thus, arbitrary profiles could be used as the 
initial conditions and the solution would develop into 
its similarity form. In our case, constant values are 
assigned to the velocity, turbulent energy and the dis- 
sipation rate. The exact values for the initial profiles 
are immaterial as long as the turbulent boundary gets 
generated. The boundary conditions used are that, 
at the free stream, the velocity reaches that of the 
fxeestream and the gradients for the turbulent energy 
and the dissipation rate are set to zero. Figures 11-14 
show the predicted velocity profile, shear stress, tur- 
bulent energy, and dissipation rate at Re$ = 1410. 
The DNS data of Spalart 11 is also shown. Again, the 
predictions from the Jones-Launder model and the 
Chien’s model are shown for comparisons. 

The computations were also carried out for larger 
Reynolds number up to Re 0 - 16000. The compar- 
ison is made between the model predictions and the 
experimental results of Wieghardt and Willmann 14 
for zero pressure gradient turbulent boundary layers. 
Figure 15 shows the results of the skin friction coef- 
ficient as a function of In figure 16, the velocity 
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Figure 3: Mean velocity profile for turbulent channel 
flow at Re r = 180. 



Figure 5: Turbulent kinetic energy for turbulent 

channel flow at Re T = 180. 



Figure 4: Turbulent shear stress for turbulent channel Figure 6: Turbulent dissipation rate for turbulent 

flow at Re r = 180. channel flow at Re r = 180. 
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Figure 7: Mean velocity profile for turbulent channel 
flow at Re r = 395. 


Figure 9: Turbulent kinetic energy for turbulent 

channel flow at Re r = 395. 
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Figure 11: Mean velocity profile for ZPG turbulent Figure 13: Turbulent kinetic energy for ZPG turbu- 
boundary layer at Reg = 1410. lent boundary layer at Re e = 1410. 




Figure 12: Turbulent shear stress for ZPG turbulent Figure 14: Turbulent dissipation rate for ZPG turbu- 
boundary layer at Reg = 1410. lent boundary layer at Reg = 1410. 
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Figure 15: Skin friction for ZPG turbulent boundary 
layer. 


Figure 16: Velocity profile for ZPG turbulent bound- 
ary layer at Re$ = 8900. 


profile at Re $ = 8900 is presented. The predicted 
value of Cf is a few percent higher than that of the 
experiment. The predicted velocity profile agrees well 
with the experiment. 

When the turbulent boundary layer is subject to 
a pressure gradient, the similarity solution ceases to 
exist. In this case, accurate description of the ini- 
tial conditions for the velocity profile and the profiles 
of turbulence quantities (turbulent energy and dissi- 
pation rate) are very important. While experiment 
could provide the velocity profile at the upstream lo- 
cation, information on the turbulence quantities (es- 
pecially the dissipation rate) is hardly available. In 
this study, the issue of the initial condition is dealt 
with in the following manner. We assume that the 
turbulent boundary layer develops under zero pres- 
sure gradient until it passes into the working section 
of the wind tunnel, where the experimental measure- 
ments are made. The connecting point between this 
virtual flat plate boundary layer and the real bound- 
ary layer with the pressure gradient is determined by 
the value of Rc$ which is found by the experiment. 

The boundary conditions are specified in the same 
way as in the case of a flat plate boundary layer. At 
the wall, both the velocity and the turbulent energy 
are equal to zero while the dissipation rate is given by 
equation (12). At the free stream, zero gradients are 


assigned to the turbulent energy and its dissipation 
rate. The mean velocity approaches that of the free 
stream, which is determined by the experiment and 
is related to the pressure gradient of the flow. 

The turbulent boundary layer studied by Herring 
and Norbury 13 was chosen as a test case for the tur- 
bulent boundary layer with a favorable pressure gra- 
dient. At the first point of the working section of the 
experiment, Re$ = 3400. Thus, profiles of the mean 
velocity and turbulent quantities (k and e) of a flat 
plate boundary layer at Re e = 3400 are used to pro- 
vide the initial conditions. With the initial conditions 
given, the calculation of the boundary layer is then 
carried out downstream. The result for the skin fric- 
tion coefficient is shown in figure 17, and the result 
for the mean velocity at x = 4ft is shown in figure 18. 
The distances, x and y, are in physical units while the 
skin friction and the velocity are normalized by the 
freestream velocity at the streamwise location under 
consideration. It is seen that predictions for both the 
velocity profile and the skin friction are excellent. 

The turbulent boundary layer studied by Samuel 
and Joubert 14 was chosen as a test case for the tur- 
bulent boundary layer with an adverse pressure gra- 
dient. In this flow, the initial conditions are specified 
at x = 0.855m, where Re# = 5470. First, the pre- 
dicted skin friction is shown in figure 19 along with 


8 





Figure 17: Skin friction for FPG turbulent boundary 
layer* 



Figure 18: Velocity profile for FPG turbulent bound- 
ary layer at X — 4.0ft. 


the experimental values. The predicted skin friction 
is higher than that from the experiment, especially 
in the latter portion of the of the experimental test 
section where the adverse pressure gradient is fairly 
strong. The predicted velocity profiles at two differ- 
ent locations ( X — 2.4t 71 and X = 3.04771, respec- 
tively.) are shown in figure 20, along with the exper- 
imental data. It is seen that the velocity profiles are 
in reasonable agreement with the experimental data. 

5. Conclusions 

By introducing the parameter R> we have success- 
ful removed the y dependence of the damping func- 
tion in YS. Now, the current model is free from the 
three deficiencies of other existing near wall k - e 
models mentioned in the introduction. First, the pro- 
posed model uses the same set of model constants as 
that used in the Standard k -€ Model and away from 
the wall the proposed model will reduce to the Stan- 
dard k - € Model. Thus, the proposed model would 
be applicable for both near wall turbulence and high 
Reynolds number turbulence. Second, the proposed 
model uses a time scale which has the Kolmogorov 
time scale as its lower bound. By using this time 
scale to reformulate the dissipation equation, the sin- 
gularity in the dissipation rate equation of the Stan- 
dard k-e Model is removed as the wall is approached 
and the equation can be integrated to the wall. This 
renders the introduction of pseudo-dissipation unnec- 
essary. Third, the proposed model uses R as its in- 
dependent variable in the damping function. This al- 
lows the model to be used in more complicated flow 
situations, flows with separation, for example. In ad- 
dition, this makes the proposed model Galilean and 
tensorial invariant. Since all the quantities in the 
proposed model are given in the local variables, the 
proposal model is very suitable for general purposed 
CFD code with unstructured grid. 

Turbulent channel flows at different Reynolds num- 
bers and turbulent boundary layers with/without 
pressure gradient are calculated using the present 
model. At low Reynolds number, the comparison be- 
tween the DNS data and the present model is found 
to be excellent. At higher Reynolds numbers, the 
velocity profiles are well predicted in all cases. How- 
ever, the predicted skin friction does not respond ad- 
equately with the pressure gradient. All the other 
existing k — e models suffer from the same deficiency, 
as pointed out by Wilcom 1 5. Currently, effort is be- 
ing made to improve the model in this aspect. 

It should be mentioned that the model is computa- 
tionally robust. Arbitrary initial profiles can be used 
for turbulent channel flows and flat plate boundary 
layers when similarity solutions exist. The predicted 
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0.0040 



Figure 19: Skin friction fox APG turbulent boundary 
layer. 



Figure 20: Velocity profiles for APG turbulent 

boundary layer at X = 2.4m and X = 3.04m. 


solution is also found to be quite insensitive to the 
number of grid points. 
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proposed model is then Galilean and tensorial invariant. The model constants used are the same as in the high 
Reynolds number Standard k-e Model. Thus, the proposed mode) will be also suitable for flows far from the wall 
Turbulent channel flows and turbulent boundary layer flows with and without pressure gradients are calculated 

T T„ d ' reCt numerical solutions and experiments show that the model predictions are 
excellent for turbulent channel flows and turbulent boundary layers with favorable pressure gradients, good for 
Uirbulent boundary layers with zero pressure gradients, and fair for turbulent boundary layer with adverse pressure 
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